#include "opt.h"

float insert32(float (*f)(float), float mu1, float mu2, float mu3, float epsilon){
    if(!(mu1 < mu2 && mu2 < mu3)){
        fprintf(stderr,"mu1 < mu2 < mu3 ???\n");
        exit(1);
    }
    if(f(mu1) < f(mu2) || f(mu3) < f(mu2)){
        fprintf(stderr, "你确定f(mu2)最小吗?\n");
        exit(1);
    }
    float Phi1 = f(mu1), Phi2 = f(mu2), Phi3 = f(mu3);
    float mu = 0, Phi, A;

    do{

        //step2
        A = 2 * (Phi1 * (mu2 - mu3) + Phi2 * (mu3 - mu1) + Phi3 * (mu1 - mu2));
        if(A == 0){
            mu = mu2; Phi = Phi2;
            break;
        }
        //step3
        mu = (Phi1 * (mu2 * mu2 - mu3 * mu3) +\
            Phi2 * (mu3 * mu3 - mu1 * mu1) + \
            Phi3 * (mu1 * mu1 - mu2 * mu2)
            ) / A;
        if(mu < mu1 || mu > mu3){
            mu = mu2; Phi = Phi2;
            break;
        }

        //step4
        Phi = f(mu);
        printf("%.3f;  %.3f;  %.3f; %.3f; %.3f; %.3f; %.3f; %.3f; %.3f\n", \
        mu1, mu2, mu3, Phi1, Phi2, Phi3, mu, Phi, fabs(mu2 - mu));

        if(fabs(mu - mu2) < epsilon)break;
        //step 5
        if(mu < mu3 && mu > mu2){
            if(Phi < Phi2) {mu1 = mu2; Phi1 = Phi2; mu2 = mu; Phi2 = Phi;}
            else {mu3 = mu; Phi3 = Phi;}
        }else{
            if(Phi < Phi2){mu3 = mu2; Phi3 = Phi2; mu2 = mu; Phi2 = Phi;}
            else{mu1 = mu; Phi1 = Phi;}
        }
    }while(1);
    printf("mu = %f, Phi = %f\n", mu, Phi);
}

float q1function(float x){
    return exp(x) - 5 * x;
}

float q2function(float x){
    return x * x * x - 2 * x * x + 1;
}

int main(int argc, char **argv){
    insert32(q1function, 1, 1.5, 2, 0.001);
    insert32(q2function, 0, 1, 3, 0.001);
}


